Calling External C functions from Mathematica

Steve Wilson
Wolfram Research, Inc.

MathLink allows the programmer to access functionality contained in libraries or other programs from the Mathematica environment.  This exercise will familiarize the programmer with the necessary techniques to write MathLink C code, template files, and integrate external programs into the Mathematica environment.  The exercises contain Answer sections where appropriate.  Please try to work out the answer to the exercise on your own before examining the provided solution. Credit goes also to Todd Gayley for some of the content in the exercises.

Component Description

For the purposes of this exercise you will use the following elements from the MathLink Developer Kit:  mathlink.h, libML32i3.a, mprep and mcc.  The header file mathlink.h contains all the data type and function definitions for the MathLink C Application Programming Interface.  Here are a few of the functions from mathlink.h:

MLDECL( MLEnvironment, MLInitialize,   ( MLParametersPointer p));
MLDECL( MLINK,         MLOpenArgcArgv,  ( MLEnvironment ep, char** argv, char** argv_end, long* errp));
MLDECL( mlapi_result,   MLPutInteger,       ( MLINK mlp, int i));
MLDECL( mlapi_result,   MLPutLongInteger,   ( MLINK mlp, long l));
MLDECL( mlapi_result,   MLGetInteger,       ( MLINK mlp, int* ip));
MLDECL( mlapi_result,   MLGetLongInteger,   ( MLINK mlp, long* lp));
.
.
.
MLDECL( void,          MLClose,           ( MLINK mlp));
MLDECL( void,          MLDeinitialize, ( MLEnvironment env));

Don't be bothered by the MLDECL() C preprocessor macro.  You can think of the macro as having the following definition:

#define MLDECL(returntype, functionname, arguments) \
	returntype functionname(arguments)

The programmer creates a MathLink program by using the MathLink API in a C/C++ program and linking that program with the MathLink library; the library of C functions that implement the MathLink protocol.

For simple MathLink programs, programmers can use the MathLink template facility to avoid having to write most of the MathLink interaction code.  Programmers simply create a template file that describes the input and output parameters to external C functions, the Mathematica calling convention for the function and allows the MathLink tool mprep to construct the necessary glue code based on the contents of the template file.

In order to compile your solutions to the exercises below, you will need to set yourself up to make use of mcc the compiler driver included in the Unix and Mac OS X MathLink Developer Kits.  To use mcc you need to open a terminal window on your machine.  (If you are using Mac OS X, you need to launch /Applications/Utilities/Terminal.app, on Linux you can just create a new xterm, and on Windows you will need to open a window that gives you access to an OS X or Linux machine.  You can complete all the exercises on a Windows machine, but that requires a compiler such as Microsoft's Visual Studio .NET.  The instructions below typically refer to the Unix/OS X environments for the sake of brevity.)  Determine the type of interactive shell you are using by running the shell command:

prompt> echo $SHELL

For /bin/tcsh and /bin/bash use the alias command to set up an mcc alias.  For OS X and tcsh:

prompt> alias mcc '/Applications/Wolfram/Mathematica\ 6.0.app/SystemFiles/Links/MathLink/DeveloperKit/MacOSX/CompilerAdditions/mcc'

For /bin/bash or /bin/sh and OS X use:

prompt> alias mcc="/Applications/Mathematica\ 6.0.app/SystemFiles/Links/MathLink/DeveloperKit/MacOSX/CompilerAdditions/mcc"

On Unix machines you can find the MathLink Developer Kit in /usr/local/Wolfram/Mathematica/6.0/SystemFiles/Links/MathLink/DeveloperKit/$SystemID.

The AddTwo Exercise

The MathLink Developer Kit contains an example of an external MathLink program called addtwo.  Addtwo implements a simple function that adds two integers and returns the result.  The Developer Kit contains the two files that make up the addtwo source code(see $TopDirectory/AddOns/MathLink/DeveloperKit/$SystemID/MathLinkExamples).

Here is the C source code for addtwo:

#include "mathlink.h"

int addtwo( int i, int j);

int addtwo( int i, int j)
{
        return i+j;
}

int main(int argc, char* argv[])
{
        return MLMain(argc, argv);
}

The main() function used here works with the Unix or Mac OS X platforms.  The following function will work for Windows:

int PASCAL WinMain( HINSTANCE hinstCurrent, HINSTANCE hinstPrevious, LPSTR lpszCmdLine, int nCmdShow)
{
        char  buff[512];
        char FAR * buff_start = buff;
        char FAR * argv[32];
        char FAR * FAR * argv_end = argv + 32;

        hinstPrevious = hinstPrevious; /* suppress warning */

        if( !MLInitializeIcon( hinstCurrent, nCmdShow)) return 1;
        MLScanString( argv, &argv_end, &lpszCmdLine, &buff_start);
        return MLMain( (int)(argv_end - argv), argv);
}

These main() functions simply call the function MLMain() that mprep creates to handle the establishment of the MathLink connection and manage the data transfer.

At this point you should note the simplicity of the addtwo function.  Addtwo is a simple function that involves no direct calls to the MathLink API or any other complicated MathLink functionality.  Often times external MathLink programs do not need any greater complexity.  Programmers can write the functions that accomplish their purpose and allow mprep to handle all the calls to the MLPut*() and MLGet*() functions(MLPutInteger(), MLPutDouble(), MLGetLongInteger(), MLGetString(), etc...).

The addtwo template file contains the template directives that give mprep instructions on how to generate the necessary glue code.

Here is the template file addtwo.tm

extern int addtwo P(( int, int));

:Begin:
:Function:       addtwo
:Pattern:        AddTwo[i_Integer, j_Integer]
:Arguments:      { i, j }
:ArgumentTypes:  { Integer, Integer }
:ReturnType:     Integer
:End:

:Evaluate: AddTwo::usage = "AddTwo[x, y] gives the sum of two machine integers x and y."

A MathLink template file may contain any number of template entries and arbitrary C code.  Each entry in the template file details the calling conventions for a function defined in the C source code.  The addtwo program contains only one function addtwo() so the template file has only one entry.  Mprep delineates template entries using the :Begin: and :End: tokens.  The :Function: token indicates the name of the C routine to pair with the template entry.  The :Pattern: token contains the function definition for the Mathematica side of the interface.  The value of the :Pattern: token effectively becomes the left hand side of the function definition much as if you did the following:

AddTwo[i_Integer,j_Integer]:=i+j

The :Arguments: token describes the arguments to the C function in the Mathematica language.  You can think of the :Arguments: token as a potential preprocessing treatment of the arguments.  Instead of just passing the arguments to AddTwo[] through to the C function addtwo(), the :Arguments: line could just as well be:

:Arguments: { i^3, 3*j }

The :ArgumentTypes: and :ReturnType: tokens take special values recognized by mprep to describe the type of input values and the return type supported by the C function addtwo.

The :ArgumentTypes: token supports the following values:

The :ReturnType: token supports the following values:

The :Evaluate: token indicates a line of arbitrary Mathematica input for the Kernel to evaluate.  Thus as in the addtwo template file a programmer can use the :Evaluate: line to create usage messages, error messages or any number of other functionality to go along with a function template.

To compile the code in the template file and the code in the addtwo source file, you can use several methods.  The easier (on Unix or Mac OS X) method makes use of the complier script mcc found in $TopDirectory/AddOns/MathLink/DeveloperKit/$SystemID/CompilerAdditions. 

Here is how to use mcc to generate the addtwo binary.

mcc -o addtwo addtwo.c addtwo.tm

The mcc script calls mprep and the system C compiler to generate addtwo.  If instead you need to use a specialized compiler or have other intermediate steps you can instead run mprep and the compiler from a Makefile.  You can find an example of a Makefile set up in $TopDirectory/SystemFiles/Links/MathLink/DeveloperKit/$SystemID/MathLinkExamples/Makefile.

The following demonstrates the mprep invocation to generate the addtwo glue C code.

mprep -o addtwo.tm.c addtwo.tm

Finally, to use the new addtwo program in the Mathematica Kernel you must install the functionality into the kernel using the Kernel function Install[].

link=Install["~/addtwo"]

LinkPatterns[link]

You can know call the function AddTwo[] and get a result.

AddTwo[10,20]

To remove the function AddTwo[] from the Kernel session and stop using the external program use Uninstall[].

Uninstall[link]

Now, use the terminal program to create a new terminal environment, create the addtwo files from the content in this notebook, compile the addtwo example program, and Install[] the result in the Mathematica session.

The AddIntegerList Exercise

Building on the previous exercise, want to now develop a simple function AddIntegerList[] that will take a list of integers as an argument, compute the sum of the elements and return the result.  Add the new function templates and C code to the addtwo files.

For this exercise you will want to specify the :ArgumentTypes: for AddIntegerList[] as:

:ArgumentTypes: { IntegerList }

With an argument of IntegerList, mprep expects the C function to take the two arguments:

void function(int *list, long n);

The argument list is a pointer to an array of int types holding the value of the list.  The argument n holds the total number of elements in list.

Now create the function template and C code resulting in the following:

AddIntegerList[{34,23,43,12}]

Answer

Additions to addtwo.tm

extern int addintegerlist(int *, long);

:Begin:
:Function: addintegerlist
:Pattern: AddIntegerList[list_List]
:Arguments: { list }
:ArgumentTypes: { IntegerList }
:ReturnType: Integer
:End:

Additions to addtwo.c

int addintegerlist(int *, long);

int addintegerlist(int *list, long n)
{
	int sum = 0, i;
	
	for(i = 0; i < n; i++){
		sum += *(list + i);
	}
	
	return sum;
}

The IntegerListTranspose Exercise

Building on the previous exercise, we want to now develop a simple function that imports a matrix of integers, transposes the matrix and returns the result. We will call the function IntegerListTranspose[]. Add the new template entries and C code to the addtwo files.

Matrices in Mathematica are simply lists of lists as we can see from the very simple matrix (({{1, 2}, {3, 4}})) in InputForm notation:

(({{1, 2}, {3, 4}}))//InputForm

Our IntegerListTranspose[] function should behave in this manner:

MatrixForm[IntegerListTranspose[(({{2, 3}, {4, 5}}))]]

MatrixForm[IntegerListTranspose[(({{6, 1, 3}, {3, 4, 9}, {8, 7, 5}}))]]

In order to transfer multi-dimensional arrays in an external mprep program we will have to make use of the MathLink API, specifically the functions MLGetIntegerArray(), MLPutIntegerArray() and MLDisownIntegerArray().  The C standard does not define matrix types as a convenient data type that can easily pass between functions and instead relies on arrays and pointers to make multidimensional collections of data.  In order to receive the multidimensional array we need to set the :ArgumentTypes: token to:

:ArgumentTypes: Manual

In order to return a multidimensional array from our function we need to set the :ReturnType: token to:

:ReturnType: Manual

These types tell mprep that the function will handle its arguments from the Kernel and its return type by using the MathLink API.  Because both the argument types and return types are Manual our new function prototype has to have following form:

void integerlisttranspose(void);

The following snippet of code illustrates how to use MLGetIntegerArray() to read the matrix from the Kernel:

long *dimensions; 
char **heads; 
long depth;
int  *data;

MLGetIntegerArray(stdlink, &data, &dimensions, &heads, &depth);
...
MLDisownIntegerArray(stdlink, data, dimensions, heads, depth);

In all mprep programs the link to the kernel can be accessed through the MLINK object stdlink.  The variable depth stores the number of dimensions of the array.  The dimensions variable has a list of depth elements where each element gives the length of the matrix in each dimension.  The heads variable points to a list of strings that represent the Head[] of each dimension.  The data variable then points to an array of integers holding the actual matrix values.

MathLink manages all the memory allocation when it gives you handles to the data, dimensions and heads from a multidimensional array.  MLGetIntegerArray() sets the data, dimensions and heads pointers up so that they access the data as retrieved and formatted by the MathLink library.  However, once you have finished using the contents of the array, you need to let the MathLink library know that it can free the allocated resources associated with the data.  MLDisownIntegerArray() effectively tells MathLink that you no longer need to use the memory that previously held the array.

The following code snippet illustrates the use of MLPutIntegerArray() to transfer an array of integer data to the Kernel.

long dimensions[2];
char *heads[2] = {"List", "List"};
long depth = 2;
int  *mat;
int  i,j, k = 0;
int  n = 3;
mat = (int*) calloc(n * n, sizeof(int));

for(i=0; i<n; i++)
	for(j=0; j<n; j++)
		mat[i + j] = k++;

dimensions[0] = dimensions[1] = n;

MLPutIntegerArray(stdlink, mat, dimensions, heads, depth);

free(mat);

Using MLPutIntegerArray() (or any of the MLPut*Array() functions) you do not need to create a heads array if all the heads of the matrix are List.  Simply pass MLPutIntegerArray() a NULL pointer and the library will add the heads for you.

At this point you have all the tools necessary to create a the new function integerlisttranspose and its template entry.

Answer

Additions to addtwo.tm

extern void integerlisttranspose(void);

:Begin:
:Function: integerlisttranspose
:Pattern: IntegerListTranspose[list_?MatrixQ]
:Arguments: { list }
:ArgumentTypes: { Manual }
:ReturnType: Manual
:End:

Additions to addtwo.c

/* Be sure to add the following to the top of addtwo.c */
#include <stdlib.h>
#include "mathlink.h"

void integerlisttranspose(void);

void integerlisttranspose(void)
{
	long *dimensions; 
	char **heads; 
	long depth;
	int  *data;
	int *transposeddata;
	long transposeddimensions[2];
	long i, j, size;
	
	MLGetIntegerArray(stdlink, &data, &dimensions, &heads, &depth);
	
	/* Calculate the size needed for the transposed array */
	size = dimensions[0] * dimensions[1];
	
	/* Create the transposed array storage */
	transposeddata = (int *)calloc(size, sizeof(int));
	
	for(i = 0; i < dimensions[0]; i++){
		for(j = 0; j < dimensions[1]; j++){
			transposeddata[i + j * dimensions[0]] = data[j + i * dimensions[1]];
		}
	}
	
	transposeddimensions[0] = dimensions[1];
	transposeddimensions[1] = dimensions[0];
	
	MLPutIntegerArray(stdlink, transposeddata, transposeddimensions, heads, 2);
	
	free(transposeddata);
	MLDisownIntegerArray(stdlink, data, dimensions, heads, depth);
}

The IntegerListTranspose Variation Exercise

For an extra challenge use Flatten[] on the Mathematica side of IntegerListTranspose[] to condense the multidimensional matrix into a simple list.  You can then use the IntegerList :ArgumentTypes: value to get the list.  You will also want to transmit the dimensions of the original matrix in the arguments to the C function.  Call the new function IntegerListTranspose2[]

Answer

Addtions to addtwo.tm

extern void integerlisttranspose2(int *, long, int *, long);

:Begin:
:Function: integerlisttranspose2
:Pattern: IntegerListTranspose2[list_?MatrixQ]
:Arguments: { Dimensions[list], Flatten[list] }
:ArgumentTypes: { IntegerList, IntegerList }
:ReturnType: Manual
:End:

Additions to addtwo.c

void integerlisttranspose2(int *, long, int *, long);

void integerlisttranspose2(int *dimensions, long dimn, int *data, long n)
{
	int *transposeddata;
	long transposeddimensions[2];
	int i, j;
	long size;
	
	size = dimensions[0] * dimensions[1];
	
	transposeddata = (int *)calloc(size, sizeof(int));
	
	for(i = 0; i < dimensions[0]; i++){
		for(j = 0; j < dimensions[1]; j++){
			transposeddata[i + j * dimensions[0]] = data[j + i * dimensions[1]];
		}
	}
	
	transposeddimensions[0] = dimensions[1];
	transposeddimensions[1] = dimensions[0];
	
	MLPutIntegerArray(stdlink, transposeddata, transposeddimensions, (char **)0, 2);
	
	free(transposeddata);
}